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Abstract 

The gradient crystal plasticity framework of Wulfinghoff et al. [53], incorporating an 
equivalent plastic strain 7eq and grain boundary yielding, is extended with grain bound¬ 
ary hardening. By comparison to averaged results from many discrete dislocation dy¬ 
namics (DDD) simulations of an aluminum type tricrystal under tensile loading, the new 
hardening parameter of the continuum model is calibrated. Although the grain bound¬ 
aries (GBs) in the discrete simulations are impenetrable, an infinite GB yield strength, 
corresponding to microhard GB conditions, is not applicable in the continuum model. A 
combination of a finite GB yield strength with an isotropic bulk Voce hardening relation 
alone also fails to model the plastic strain profiles obtained by DDD. Instead, a finite GB 
yield strength in combination with GB hardening depending on the equivalent plastic 
strain at the GBs is shown to give a better agreement to DDD results. The differences 
in the plastic strain profiles obtained in DDD simulations by using different orientations 
of the central grain could not be captured. This indicates that the misorientation depen¬ 
dent elastic interaction of dislocations reaching over the GBs should be included in the 
continuum model, too. 

Keywords: Grain boundary hardening, grain boundary yielding, strain gradient 
plasticity, discrete dislocation dynamics. 
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Microstructural characteristics influence the material behavior of multicrystalline metals. 
Their mechanical response is controlled by the dislocation microstructure, their interac¬ 
tions, and their complex interplay with other microstructural characteristics like, e.g., 
precipitations and grain boundaries (GBs) [9]. Dislocations can, e.g., pile up at GBs, be 
absorbed into the GBs, or be transmitted into the adjacent grains [38]. To account for 
these GB mechanisms in continuum models is still an open challenge, and only a limited 
number of three-dimensional models exist (e.g., [27]). In the present work, the focus is on 
GB effects on the local material behavior. On that account, a phenomenological model 
for GB hardening is incorporated into a gradient plasticity (GP) framework. Regarding 
its prediction of the plastic strain close to GBs, it is evaluated by comparison to discrete 
dislocation dynamics (DDD) results. 

Gradient plasticity theories are commonly used to take into account microstructural char¬ 
acteristics in continuum modeling. In order to model, e.g., size effects [16], an internal 
length scale is introduced. Therefore, an additional (defect energy) contribution to the 
Helmholtz free energy of the bulk is usually considered. This contribution is formu¬ 
lated, for example, using the gradients of plastic slip or the plastic deformation gradient 
(e.g., [25]). Such an extension requires equilibrium conditions for the additionally aris¬ 
ing microstresses. Thus, higher order boundary conditions (BCs) are needed, as well. 
Choices for these can, but not necessarily must, result from the requirement of achieving 
a null-expenditure of power by the microstresses. In [31], two such choices are discussed. 
These give upper and lower bounds for the defect transfer restrictions imposed by the 
boundary / grain boundary tractions. They are referred to as “microscopically free” 
(microfree) and “microscopically hard” (microhard) conditions. The microfree case cor¬ 
responds to vanishing tractions for the microstresses, i.e., the tractions do not impede the 
flow of dislocations which itself is modeled, e.g., by plastic strain rates. An unrestricted 
flow of dislocations across the boundary corresponds to continuous plastic strain rates. 
Whether these are continuous also depends on the material model applied and the me¬ 
chanical boundary conditions. For the opposing microhard condition, the slip rates are 
set to vanish on the boundary. This condition corresponds to a non-passing restriction on 
dislocations. Both BGs can also be applied for the case of only one accumulated plastic 
slip. In [30], this is discussed, e.g., for the gradient theory of Aifantis [2, 3]. However, 
these two bounds for BGs are not applicable to all types of GB behavior. Van Beers et al. 
[46], e.g., give an overview about the possible interactions between dislocations and GBs. 
One approach to achieve GB behavior inbetween the above described “simple” conditions 
is the introduction of a GB energy (see, e.g., the quadratic energies in [6, 24], and the 
more general forms in [22, 48]). Grain boundary hardening can be modeled with such 
energy approaches (see also [5, 47]). 

The consideration of GBs in gradient plasticity models requires to specify the continuity 
of quantities like plastic strain and microstresses at these interfaces. In Aifantis et al. 
[5], continuity of the plastic strain is assumed. Jumps in its gradient, however, are al¬ 
lowed for. This model has shown good agreement with experiments for a one-dimensional 
case. In Gurtin and Needleman [31], jumps in the Burgers vector flow across the GB 
are considered. The higher-order stresses are, however, assumed continuous across the 
GB. For a detailed overview on modeling approaches with regard to the continuity of 
higher-order quantities see [47]. In the field theory of defects by Fressengeas et al. [23], 
the tangential continuity of different elastic/plastic tensors across an interface is derived 
from the conservation of Burgers and Frank vectors, and compatibility conditions are 
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established at multiple junctions. 

Naturally, continuum models and the employed hardening relations are intended to be 
calibrated to experiments. Latent hardening models (e.g., [21, 34, 36]) are commonly 
used (e.g., [7, 8]). They take hardening occurring in other (than primary) slip systems 
due to slip in primary slip systems into account. The Mises-Hill framework of Gurtin and 
Reddy [32] considers self hardening (i.e., slip systems harden due to their own plastic slip) 
and latent hardening. A hardening rule is introduced on the individual slip systems and 
takes accumulated plastic slips into account. Furthermore, hardening approaches based 
on scalar plastic strain measures have been developed (see [13] for a brief overview). Such 
approaches allow for rapid parameter calibration with a manageable number of experi¬ 
ments. 

However, it is challenging to develop appropriate experimental settings that can be easily 
varied to investigate several test cases. To overcome these experimental challenges, DDD 
simulations are used in the present work. The physically detailed modeling of plastic¬ 
ity by DDD results in information about the effects of the dislocation interactions. In 
such approaches (e.g., [15, 33, 39, 41, 43, 45]), the plastic response is directly computed 
from individual dislocation motions. Dislocations are described as elastically interacting 
line defects, and are discretized as polygons. Physical mechanisms of dislocation glide, 
cross-slip and reactions are treated with constitutive rules. In the DDD framework of the 
present work (Weygand et al. [50, 51], Siska et al. [49]), GBs confine dislocation motion 
to remain within the respective grains. Nevertheless, total transparency of the GBs to 
stress and displacement fields of each dislocation (and, thus, the elastic interaction of dis¬ 
locations across GBs) is preserved. From the results of DDD simulations, conclusions can 
be drawn to incorporate the dislocation behavior in the GP model phenomenologically. 
Bardella et al. [11], e.g., use a strain gradient extension of conventional latent hardening 
and benchmark it with DDD results for a simple shear problem. In Aifantis et al. [4], 
a tricrystal tensile setting of microsize is used to compare results of a GP model with 
DDD results. The GB modeling in the DDD simulations is the same as in the present 
work, i.e., dislocation glide is restricted to remain within the grains. Microhard bound¬ 
ary conditions are utilized in the used GP model of [4]. However, in a reevaluation of 
the data and GP theory of [4], Zhang et al. [56] showed that a much better agreement 
between DDD and GP results is obtained by using GB yielding in the GP model. 

The theoretical basis for the GP model of the work at hand is the framework by Wulf- 
inghoff et al. [53] which supplements [54] (see [55] for the time integration algorithm) 
by a GB energy and a GB yield condition. In the energetic framework used, the GB 
energy models the storage of energy due to the accumulation of defects at the GBs. The 
GB yield strength models the resistance of the GBs against plastic flow, necessary to 
match the strain profiles obtained by DDD in this work. For numerically efficient com¬ 
putations, the gradient contribution to the free energy is formulated with respect to a 
scalar equivalent plastic strain 7 eq, instead of considering gradients of all plastic slips 
individually. Thus, a pragmatic engineering approach rather than a physical approach 
is taken. Usually, the magnitude of size effects modeled with GP is associated to the 
internal length scale used in models. Such a length scale is introduced in the mentioned 
framework via a quadratic defect energy (see [53], and Section 2.2.3 for details). Increas¬ 
ing the GB yield strength in this model, however, significantly intensifies the magnitude 
of size effects observed. Ghanging the internal length scale has a significantly smaller 
effect, see [53]. Additionally, in [53], an isotropic Voce hardening law for the bulk, also 
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posed in terms of 7 eq, is used. This phenomenologically accounts for hardening associ¬ 
ated to mechanisms like dislocation trapping. The combination of these approaches has 
shown to give good agreement with the size effect results from a tensile test experiment 
on microwires [53]. In addition, it has been shown in [57] that the 7 eq-distributions, 
obtained from simulations of single-crystalline gold microwires, are compatible with the 
experimentally determined kernel average misorientation distributions in cross sections 
of these microwires. 

Due to numerical reasons, a micromorphic approach, suggested by Forest [18], is used 
for the implementation of the theory at hand. In the micromorphic approach, additional 
variables are introduced as additional internal degrees of freedom. Micromorphic the¬ 
ories belong to the class of generalized continuum theories. It can be formally shown 
that gradient theories, which also belong to this class, are special cases of micromorphic 
theories (see [20], and also [37] for an overview on dislocations in generalized continua 
frameworks). In [18], a general framework for the micromorphic approach is outlined. 
It includes balance equations governing the micromorphic degrees of freedom, boundary 
conditions, and higher-order stresses (see also the references in [18]). An overview of 
applications to, e.g., elasticity and gradient plasticity is also given. In this context, it is 
described in detail that models formulated with the micromorphic approach are related 
to existing gradient models. Therefore, the micromorphic variable is constrained to be 
equal to its (macro) counterpart. This constraint can be imposed by a penalty term in 
the free energy. The constrained micromorphic approach yields models that, e.g., belong 
to the class of gradient of internal variable models [40]. For instance, the gradient theory 
by Gurtin [26] can be regarded as a constrained micromorphic theory [18]. 

In the work at hand, the GP approach by [53], implemented within a constrained mi¬ 
cromorphic approach [18], is extended by GB hardening. Fully threedimensional GP 
simulations are performed. The GP results are compared to DDD results. This compar¬ 
ison includes GP results obtained by using Voce hardening on one hand, and GP results 
obtained by using GB hardening on the other hand. Different crystal orientations and, 
consequently, different slip system mismatches between adjacent grains, are considered. 

Notation. A direct tensor notation is preferred throughout the text. Vectors and 2nd- 
order tensors are denoted by bold letters, e. g., by a or A. A linear mapping of 2nd-order 
tensors by a 4th-order tensor is written as A = C[S]. The scalar product and the dyadic 
product are denoted, e.g., hy A ■ B and A® B, respectively. The composition of two 
2nd-order tensors is formulated by AB. Matrices are denoted by a hat, e. g., by e. 

2. Gradient Plasticity Model 

2.1. Motivation for the Grain Boundary Hardening Extension 

The presence of GBs in microstructured materials leads to dislocation pile-ups that in¬ 
fluence the overall work hardening. In many crystal plasticity continuum models, the 
dislocation-induced hardening mechanisms are modeled by, e.g., isotropic hardening rela¬ 
tions for the bulk material. When solely the overall mechanical properties are of interest, 
this procedure is able to achieve good results. In the work at hand, however, the dis¬ 
tribution of plastic strain is evaluated, additionally. It is investigated, if the localization 
of plastic strain, resulting from pile-ups of dislocations observed at impenetrable GBs 
in DDD simulations, can be accounted for by a GP model with GB yielding. A direct 
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translation of the DDD grain boundary conditions to the GP model is not possible due 
to the coarsening in the continuum approach. In order to account for work hardening, 
at first, a bulk Voce law is used since it has given promising results in combination with 
GB yielding in a previous work [53]. It is shown here, however, that this approach is 
not able to capture the evolution of the plastic strain in the vicinity of the GBs. As an 
alternative to this bulk hardening model, a GB hardening relation is used, subsequently. 
This relation takes the equivalent plastic strain at the GBs into account. The explicit 
association of hardening to the GBs is motivated by observations in the DDD simula¬ 
tions. Dislocations entangle localized as pile-ups close to the GBs. In the context of GB 
yielding, it is also noteworthy that in Aifantis et al. [5], GB yield stresses are estimated 
based on indentation studies. The authors attributed the observed increase in hardness 
near the GB (and, thus, an increase in resolved shear stress) to dislocation pile-ups [44]. 
Keeping this in mind, the explicit consideration of the GBs in the hardening relation of 
the GP model is justihable. 

2.2. Mathematical Framework for Gradient Plasticity with Grain Boundary Hardening 
2.2.1. Basic Assumptions 

In a geometrically linear framework for deformations, the strain tensor reads e = sym(VM) 
The displacement gradient is given by Vit = d^jUiSi ® e^, in terms of the basis vec¬ 
tors {ei, 62 , 63 } in a Gartesian coordinate system. Furthermore, the plastic strain tensor 

eP = 5 ]A„M:. ( 1 ) 

Ct 

The symmetric part of the Schmid tensor reads AT® = sym(d„ 0 n„) with the slip sys¬ 
tems a = 1,2,..., N. Slip directions are indicated by da and the slip plane normals 
by tIq . It should be noted that in the following plastic slip parameters Xa are used. They 
increase monotonously, i.e., Aq > 0 in this work. Since the focus is on face-centered cubic 
(FGG) crystals, the number of slip parameters is N = 24. Gonsidering 12 slip systems, 
each with two directionally dependent slips, would be equivalent. Additive decomposi¬ 
tion of the strain tensor leads to the elastic strain tensor e® = e — e^, where e is the total 
strain tensor. The equivalent plastic strain measure is introduced as (cf. [54]) 

7 eq(A) = ^ J Aadt = ^Aa. ( 2 ) 


2.2.2. Principle of Virtual Power 

The field equations can be derived from the principle of virtual power. It states that 
the virtual power of the internal forces 5Pint equals the virtual power of the external 
forces SVext- Body forces are not present in the following derivations. The internal 
power density of the bulk, pvob is assumed to be given by 

Pvoi = o-• £-f ttC- l-| • VC- (3) 

This statement is an extension of the classic power of internal forces. It takes the power 
expended by generalized stresses, via the rate of a micromorphic field variable C and via 
the rate of its gradient, into account [19]. The stresses cr, tt, and | are work conjugate 
to £, f, and V(^, respectively, where the scalar micromorphic variable is a different 
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quantity than the magnitude ^ of the vectorial microstress The body B has an external 
boundary dB and internal boundaries denoted by the union T of all GBs. Then, the 
internal power Virit is assumed to be given in terms of the volume integral over the 
internal bulk power density Pvoi and the area integral over the power density on the 
union T of all GBs, i.e., 

'Pint= J ^voiduT J SrCda. (4) 

B r 

The GB microtraction Sp is associated to ( on the GBs. In fact, this microtraction 
is imposed in form of a jump condition for the microstresses on the GBs, as shown 
below. Furthermore, the micromorphic counterpart C of the equivalent plastic strain 7 eq 
is assumed to be continuous across the GBs, cf. [53]. Consequently, possible jumps in ( 
are neglected in the following (see also the discussion in Section 6.5). 

The external power Ve^t is assumed to consist of the following two contributions 


i^ext = J (t -11 + da. (5) 

dB 

This takes into account the power expended by tractions i and microtractions S, respec¬ 
tively, at the external boundary dB. Due to the simplification of considering an equivalent 
plastic strain 7 eq, and an additional micromorphic field variable C, with its gradient VC, 
only one microtraction S associated to the rate C is accounted for on the external bound¬ 
ary. This microtraction is imposed as a boundary condition for the microstress, as it is 
shown below. In models that account for the individual slips or dislocation densities as 
field variables it is, in principle, possible to prescribe microtractions for the individual 
slip systems. Whether this is applicable depends on the model and the employed higher- 
order quantities (e.g., [28, 31]). In the following, u = 5u and C = where 5u and 5^ 
are arbitrary virtual rates that vanish at the Dirichlet boundaries dBy^ for given {it, C}. 
The virtual power of the external forces is then given by 

dPext = I i-Su da + j SdCda, (6) 

dBt ass 

with dB = dBt U dB^. U dBu and dBt U dB^ H dBu = 0 . For the same choice of virtual 
rates, the virtual power of the internal forces, consequently, reads 

SPint = J ^(T • (5e-I- ttJC + C ■ V(5C^ du-I- J ErS(da. (7) 

B r 

Here, it should be noted that the variables and C are, a priori, chosen to be indepen¬ 
dent, i.e., = 0 (cf. [52]). With SVtnt = 5Vex_t and a standard procedure (e.g., [29]), 

it is possible to derive the field equations, the Neumann boundary conditions, and the 
GB conditions (see also [52, 54]). In Box 1, the resulting field equations are summa¬ 
rized. The classic balance of linear momentum is supplemented by a microforce balance. 
This takes into account the microstresses occurring due to the introduction of additional 
contributions to the internal power density. These terms consider the rates of the mi¬ 
cromorphic field variable and its gradient V^. Gonsequently, additional microtraction 
conditions for the gradient stress ^ at the GBs and at the external boundary supplement 
the Neumann BG of the Gauchy stress. 
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Box 1: Field equations and boundary conditions. The jump of ^ is denoted 
by 14] = 4^ ~ • The grain boundary normal points from to “+”. 


Linear momentum balance 
Microforce balance 
GB microtraction 
Neumann BCs for: Cauchy stress 
Gradient stress 


0 

= div (cr) 

VxgB 

TT 

= div(4) 

yxGB\r 

Sr 

= K1 • « 

Wx gT 

an 

= i 

on dBt 

1 n 

_ S' 

on dBs 


2.2.3. Constitutive Assumptions 

Free energy density. Following the model of [53], the free energy density in the bulk 
is assumed to have an elastic (IFe), an isotropic hardening (VFh), and a defect (Wg) 
contribution, i.e., 

1F(£, a, C, VC) = fFe(£, £P(A)) + IFh(C) + ITg(VC) + VF;,(C - 7eq(A)). (8) 

The need for additional contributions to the free energy arises from the coarsening in 
the continuum modeling of the elastic energy. However, distinct dislocation phenomena 
are modeled here using the two contributions Wh(C) and Wg(VC): on one hand, the 
isotropic hardening resulting from statistically stored dislocations, and, on the other 
hand, the influence of geometrically necessary dislocations. 

The coupling of C and 7 eq is achieved by using — 7 eq(A)) = — 7eq)^/2. This 

penalty energy ensures that C ~ 7eq for a sufficiently large penalty parameter Thus, 
in Eq. (8), the gradient extension is also performed in terms of the micromorphic field 
variable C, instead of 7 eq. The dependence of the free energy on the accumulated plastic 
slip represents a rather simple phenomenological approach. This slip does not represent a 
state variable describing the internal defect structure locally. The elastic energy density 
reads We(e,eP(A)) = (e — eP(A)) • C[£ — £P(A)]/2. Here, C denotes the elastic stiffness 
tensor. Additionally, an isotropic, rate-independent. Voce hardening relation is used via 

W4(C) = (t-<S - - Tof exp , (9) 

with the initial yield stress Tq ^ the saturation stress and the initial hardening mod¬ 
ulus 0. This hardening energy is obtained by an integration of a chosen relation for 
the hardening stress of all slip systems. Instead, it is possible to include a dissipative 
hardening contribution on each slip system. This approach would, however, lead to the 
same material behavior. The defect energy reads 

IFg(VC) = ^i^cVC • VC, (10) 

and introduces a length scale into the model by means of the defect parameter Kq (which 
is assumed to be a constant, here). The quadratic formulation of the defect energy gives 
a linear dependence of the microstress 4 on the gradient VC. This gradient could be 
interpreted as an approximative measure for geometrically necessary dislocation (GND) 
densities [53] . This simplification of representing dislocations using an equivalent measure 
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is only possible for monotonic loading processes since C ^ 0- By the quadratic choice in 


Eq. (10), the well-posedness of the boundary problem, with respect to the defect energy, 
is ensured. In this context, see the discussion of a quadratic defect energy function in [42]. 
There, however, the discussed defect energy depends on the sum over the slip gradients 
rather than on the gradient of the micromorphic counterpart of the sum over the slip 
parameters. Additionally, the gradient stresses resulting from Eq. (10) could also be 
introduced within a dissipative framework rather than in an energetic one (see, e.g., [17] 
for a fundamental discussion on the matter of dissipative and energetic frameworks). 

In addition to the bulk free energy density W, an energy density per unit surface has 
been introduced on the GBs in [53] 


B^r(C) = 5?C- 


( 11 ) 


In combination with a GB yield criterion, this leads to an explicit consideration of GB 
yielding in the continuum model, cf. also the brief discussion of GB energies in the 
introduction. 

Hardening extension of energy density on the grain boundary. The GP model is extended 
to account for hardening considering the plastic deformation of the GBs. Therefore, the 
following extension of the GB energy, Eq. (11), is proposed 


IE^(C)=sCc+iKHC^ 


( 12 ) 


where ATh is a GB hardening parameter. Contrary to Eq. (11), this GB energy leads to 
an increasing GB yield strength in the additional GB yield criterion used in this work. 
Consequently, the GB yield strength increases in dependence of f. This mechanism 
phenomenologically accounts in the continuum model for the strengthening effects caused 
by the pile-ups of dislocations at the impenetrable GBs in the discrete model. 

Dissipation. Upon neglecting thermal effects, the total dissipation reads 



(13) 


B 


r 


After exploiting T^ext = T’int, it can be summarized as 



(14) 


B 


r 


Substitution of Eq. (8) in Eq. (13) gives a bulk dissipation in Eq. (14) that reads 


P = (cr-5^We)-e-9,pWe-eP + (7r-5cIUh-5cIUx)C-5^„,IUx7eq + (|-9vcW"g)-VC > 0. 

(15) 


Introducing the abbreviation p = = —9fIU^, assuming the stresses cr, tt, and ^ 

to be purely energetic, and substituting devWe = —f leads to a reduced bulk dissipation 
inequality of the form 


V = cr ■ e’P -P7eq > 0. 
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( 16 ) 


Resulting from the micromorphic approach, the (reduced) dissipation inequality includes 
a stress p that is associated to the difference 7 eq — C (see also [18]). It is assumed that 
the bulk dissipation is induced by the dissipative shear stresses of the individual slip 
systems, e.g., [14], 

V = Y,<^o.. (17) 

Ol 

Consequently, after combining Eqs. (1) and (2) with Eqs. (16) and (17), using the micro¬ 
force balance equation from Box 1 as well as the abbreviation /3 = dcW-^, the dissipative 
shear stresses read 

Tt = 'ra + div (I) - fi, (18) 

with the resolved shear stresses = cr ■ M^, and the hardening microstress /3. 
Furthermore, the GB dissipation from Eq. (13) and Eq. (14), respectively, reads 

Vr = (Sr - Sf) C = S^ C > 0, (19) 

where Sp is the dissipative GB microtraction. The energetic GB microtraction Sp reads 


Sp = dc^Wr = Sg -I- RThC- 


( 20 ) 


Flow rules for the bulk and the grain boundary. The flow rule for the bulk is assumed to 
be of overstress type, formulated in the slip parameter rates A^, 


Aa = 7o 



7o 


Ta + div (^) - (Tq^ -h P) Y 

-- j 


( 21 ) 


with the reference shear rate %, the drag stress r^, and the rate sensitivity exponent p. 
The GB yield function is introduced as 

/r = Sp = III • n - Sp = [I] • n - (22) 


where the GB itself is assumed to deform dissipation-free. However, in principle, it 
is possible to additionally consider GB dissipation. Gontrary to a GB yield strength 
resulting from Eq. (11), the GB yield strength stemming from Eq. (12) combines a 
constant initial yield strength with an additional contribution linearly increasing with C, 
i.e., with proceeding plastic deformation. For the GB, the Kuhn-Tucker conditions read 
/r < 0, ^ > 0, C/r = 0, where the GB is assumed to behave rate-independent. 


3. Discrete Dislocation Dynamics Simulations 
3.1. Introductory Notes 

As a data basis for the parameter calibration of the GP model, a multitude of DDD 
simulations is carried out. The discrete dislocation results are obtained with the DDD 
code described in detail in [49, 50, 51], including the analysis of the elastic interactions 
of dislocations as well as the resulting plastic deformation of a sample. The framework 
allows simulating polycrystalline aggregates. For brevity, only a short description of the 
framework is given. Material parameters, interface and boundary conditions, and the 
simulation setup are outlined. The necessary ensemble averaging procedure is described 
by which simulation results directly comparable to the (non-scattering) GP simulation 
results are obtained. 
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3.2. Material Model 

In the DDD code, an individual crystallographic orientation is assigned to each grain and 
isotropic linear elasticity is used. An FCC crystal system with elastic constants mim¬ 
icking aluminum (shear modulus G = 27 GPa, Poisson’s ratio ly = 0.347) is used. The 
initial dislocation configurations are chosen such that the system is placed in the mul¬ 
tiplication controlled regime [35]. Frank-Read (FR) sources are, therefore, distributed 
randomly with respect to position and orientation. This procedure is performed un¬ 
der the restriction that each slip system in each grain contains the same number of 
sources. The initial source length is chosen to vary between 0.16 — 0.27 p-m in order 
to reduce artificiality that would be introduced into the model by the choice of a uni¬ 
form FR source length. Using this procedure, the initial dislocation density is about 
p « 7.5 X 10^^ m“^. The sources are approximately two times larger than the mean 
dislocation spacing I = l/i/p « 0.12 pm, which places the system behavior in the multi¬ 
plication controlled plasticity regime [35, 56]. Consequently, the work hardening is influ¬ 
enced by dislocation reactions, as opposed to single source controlled plasticity. Thus, a 
comparison of DDD results with results from a continuum model is feasible. 

3.3. Simulation Setup, Geometry and Boundary Conditions 

Strain rate controlled tensile tests are simulated along the (lOO)-axis (which is the x-axis 
of the in-lab frame) of a tricrystal geometry (cf. Fig. 1) with a strain rate of e = 5000 
during the simulation time t. This rate is chosen due to the small time scale at which 
DDD operates. It is assumed that the plastic behaviour is independent of the strain 
rate [43]. In order to investigate different dislocation interaction behavior across the 
GBs in the model, the crystal orientation of the central grain is rotated by an angle tp 
around the x-axis (cf. also Section 5.1). The size of the cubic grains is set to 0.75 pm 
and the BCs applied are as follows. On the two bounding planes (in loading direction) 
the displacements Ux(x = 0) =0 and Ux(x = Xmax) = x^a^et, respectively, are prescribed 
while the displacements orthogonal to the loading axis are set to zero (except for one 
special case, cf. Section 5.3). All other boundaries are traction free and dislocations are 
allowed to leave the volume, there. 

3.4.. Averaging Procedure 

Due to the discrete nature of DDD, the scatter in the results depends on the dislocation 
structure and density. The higher the density, the closer DDD results from individual 
simulations come to a continuum-like profile. However, the size of the setups considered 
in the present work is within the size effect regime, where individual dislocations control 
plasticity. Thus, a suitable averaging procedure is necessary as a prerequisite to calibrate 
GP model parameters to these results. Therefore, averaging over an ensemble of simula¬ 
tions is performed. Once all DDD simulations are carried out, the plastic strain of each 
realization is directly evaluated from the swept areas of all discrete dislocations. This 
postprocessing procedure considers slices perpendicular to the tensile axis. Subsequently, 
the DDD results from a number of simulations M are averaged, for one set of BCs and 
crystal orientations: 

1. The contribution of a slip system a to the plastic strain tensor is calculated in 
each slice of the volume V via e]) = hAal{2V) {da G -I- 0 da), where b is the 

length of the Burgers vector, and Aa the swept area of dislocations of slip system a. 
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2. The slip system contributions to the plastic strain tensor eP of a slice are su¬ 
perposed. 

3. The mean plastic strain per slice is obtained by averaging twice: 

(a) Averaging the plastic strain arithmetically, i.e., over ~ 5 DDD simulations. 

(b) Averaging the plastic strain over equivalent volumes regarding the crystallo¬ 
graphic mirror-symmetry in cc-direction at a; = O.Sxmax- Thus, the data from 
DDD is doubled and, therefore, further smoothed. This is, formally, similar 
to the averaging over all symmetry-equivalent components of the dislocation 
density tensor and symmetry-equivalent positions along the x-axis in [4]. 

The number of slices along the tensile axis is chosen to be 150, resulting in a slice thickness 
of 15 nm. This is based on the evaluation of the mean dislocation spacing in the pile- 
up close to the GBs of approximately 1/y^ « 25 nm. With 150 slices, the evaluation 
resolution is almost two times higher than the mean dislocation spacing in the pile-up, 
abundant to capture the strain gradients. The GB itself does not exhibit plastic strain, 
but the slices adjacent to it do. The averaging procedure renders the spatial distribution 
of plastic strain to be one-dimensional. All dislocations add to the produced plastic slip 
- and thus to the plastic strain - regardless of their locations in the cross sections. The 
averaged DDD results are used to calibrate the GP model parameters. 

4. Gradient Plasticity Simulation Setup 

^.1. Geometry, Boundary Conditions and Crystal Orientations 

In this section, the setup of the finite element (FE) simulations for the GP model with 
an in-house code is discussed (for details on the FE implementation see [53, 55]). In all 
simulations, a tricrystal composed of 0.75 p.m-wide cubic grains is considered under ten¬ 
sile loading with Dirichlet conditions {/S.Ux = 0.005Loi Lq = 2.25 pm, cf. Fig. 1). Lateral 
contraction on the boundary planes at x = 0, and at x = Xmaxi respectively, is prohibited 
(except for one special case, cf. Section 5.3). The restriction of lateral contraction (see 
also Fig. 1) is abbreviated by NLC in the following. At the beginning of the simulation, 
the equivalent plastic strain 7 eq (and its micromorphic counterpart C) are set to zero 
everywhere. Finite element nodes on the GB planes as well as on the boundary planes 
are assigned a GB yield strength and are set to microhard behavior, at the beginning. 
Once the (GB) yield condition at these nodes is fulfilled, they are allowed to yield. The 
motivation for assigning a GB yield strength to the boundary planes, too, is drawn from 
the BCs in the DDD simulations that lead to pile-ups of dislocations at the boundaries, 
as well (cf. Section 3.2). For details on the active set search for the (grain) boundary 
nodes in the GP simulations it is referred to [53]. 

An elastically isotropic, but plastically anisotropic aluminum-like material is considered 
throughout the following GP simulations. At first, the (lOO)-orientation of the three 
grains with FGG crystal structure is oriented aligned with the x-axis of the cartesian 
X, y, z-system. The central grain is rotated by an angle ip, subsequently (see also Sec¬ 
tion 3.3). As a first approach to investigate the interaction behavior across GBs with 
both models, three representative cases, p G {0°,5°,35°}, are selected. Thus, the ideal 
case of vanishing mismatch between the slip systems of adjacent grains (and, conse¬ 
quently, unrestricted interaction of dislocations across the GBs) is supplemented with a 
case of small mismatch (yet still strong interaction) and a case of large mismatch (weak 
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Figure 1: Equivalent plastic strain field obtained in a gradient plasticity FEM simulation 
without lateral contraction at the boundaries and finite values of the GB yield strength 
(left, depicted are the Gauss point subvolumes). DDD simulation with impenetrable GBs 
(right, dashed lines indicate the GBs). 


interaction) between adjacent grains. 

The chosen finite element mesh for the GP simulations consists of 12 x 12 x 12 elements 
for each grain, i.e., in total 25012 DOF. It is chosen as a compromise between com¬ 
putational time and accuracy. Compared to the chosen discretization, a refinement of 
the mesh with approximately twice the number of DOF yields a relative error in the 
stress-strain response at the final time-step of less than 0.01. 

^.2. Model Parameters 

The GP model parameters are calibrated such that the overall mechanical response 
matches the averaged stress-strain results of the DDD simulations. This consideration is 
not sufficient alone since distributions of the plastic strain, e.g., along a line segment of 
the whole volume, are not necessarily predicted correctly. Therefore, the local distribu¬ 
tion of plastic strain is taken into account in the calibration, too. While, in principle, the 
determination of the GP model parameters in the fitting procedure is not unique, there 
are some guidelines outlined in this work that should help in their calibration. They are 
discussed in the following with regard to the DDD results as data source. For each case 
(cf. Table 1) 

1. A least-squares fit (LSF) of the DDD stress-strain curves is obtained. 

2. The cross-section averaged plastic strain profiles along the loading axis are obtained. 
This is performed in order to ensure comparability with the (averaged) plastic strain 
profiles from DDD. The occurring differences within the cross-section distributions 
of GP results are relatively small and not as pronounced as they are in the DDD 
simulations. 

3. The Young’s modulus is calibrated in the GP model in order to match the elastic 
stiffness obtained from the LSF. 

4. The initial yield stress Tq of the GP bulk model and the initial yield strength Sq of 
the GBs as well as the initial yield strength Sq gg of the boundary (planes at x = 0 
and X = Xmax) are calibrated using the DDD-LSF and the averaged plastic strain 
profiles. Therefore, the plastic strain profiles along the loading axis of GP results and 
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DDD results are compared. On that account, plastic strain profiles are obtained at 
three, representatively chosen, fixed overall plastic strain values in the well established 
plastic regime are used, i.e., Sp S {0.001,0.002,0.003}. 

a) In case of Voce hardening: The initial hardening modulus 0 and the saturation 
stress are adjusted to the work hardening behavior of the DDD-LSF. 

b) In case of GB hardening: The GB hardening parameter is adjusted to give 
good agreement with the work hardening behavior of the DDD-LSF. Additionally, 
the evolution of the plastic strain at the GBs is taken into account, i.e., the GP 
plastic strain profiles are compared to the (averaged) DDD plastic strain profiles 
at all three overall plastic strain values Ep. 

The calibration yields a uniform Young’s modulus of A = 65 GPa for all GP simulations 
(except for one special case, cf. Section 5.3). The slightly lower value of the Young’s 
modulus, compared to the value of 72.7 GPa used in the DDD model, is due to the 
anelastic behavior (bow-out of dislocations) from the very beginning of the loading. The 
Poisson’s ratio is kept identical to the DDD simulations {v = 0.347). 

It is remarked that the defect energy parameter Kq, although in principle introducing an 
internal length scale into the model, mainly controls the elastic-plastic transition behavior 
if the GP model with GB yielding is used. Consequently, ATq = 84 x I0“® N is chosen 
such that the GP simulations show similar stress-strain results as the DDD simulations 
in this regime. Hardening resulting from the defect energy is negligible compared to the 
additional hardening relations investigated in the work at hand (cf. Fig. A.l). For all 
simulations, a reference shear rate of 70 = I 0 “^/s, a rate sensitivity exponent of p = 20 , 
and a drag stress of = 1 MPa are used. The used penalty parameter is = 10® MPa. 
All cases considered are summarized in Table 1 (see Section 5.1 for remaining Voce 
parameters). 

Table 1: Setups and model parameters of GP simulations for comparison to DDD results. 
The abbreviations NLC and LC indicate if lateral contraction is prevented or allowed 
for on the boundary planes at x = 0 and x = Xmax- In the special case LCOE, the two 
bounding grains are purely elastic. The gradient hardening contribution is negligible in 
the investigated NLC cases, see Fig. A.l. 


Name 

Angle ip 

Hardening 

"C 

“0,GB 

[N/m] 

“0,9B 

[N/m] 

ATh 

10® [N/m] 

[MPa] 

NLC35V 

35° 

Voce 

3.5 

25 

- 

30.0 

NLC35G 

35° 

GB 

1.5 

25 

1.8 

33.5 

NLC5G 

5° 

GB 

1.5 

25 

1.8 

33.5 

LCOE 

0 ° 

Gradient 

- 

- 

- 

44.0 


5. Comparison of Gradient Plasticity Results to Discrete Dislocation Dynam¬ 
ics Results 

5.1. Gradient Plasticity with Voce Hardening 

At first, the framework of [53] is used (without GB hardening) taking into account 
the calibration guidelines described above. Consequently, the averaged plastic strain 
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profiles are evaluated at three constant overall plastic strains (cf. Fig. 2, left) for the 
GP simulations. They are compared to the averaged plastic strain profiles (Fig. 2, right) 
obtained from the DDD simulations (indicated schematically by blue error bars in Fig. 2, 
left). The model parameters for the case at hand (NLC35V) can be found in Table 1. In 
addition, the saturation stress reads = 108.51 MPa, and the initial hardening modulus 
is 0 = 1075 MPa. It can be seen in Fig. 2 that, although the plastic strain profiles 
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Figure 2: NLC35V/NLC35G: Stress-strain curves of DDD and GP simulations (left). 
Distribution of cross-section averaged plastic strain along loading direction of DDD and 
GP simulations (right), obtained at fixed overall plastic strain values (left). Boundary 
conditions with restricted lateral contraction, central grain rotated 35° around loading 
axis. 


for ep = 0.001 are in good agreement as a consequence of the calibration of the GB yield 
strength, the subsequent evolution can not be accounted for by the GP simulations with 
Voce hardening. At the GBs, significant deviations occur which are caused by an obvious 
limitation of this approach to account for the observed accumulation of dislocations at 
the GBs in the discrete simulations. 

5.2. Gradient Model with Grain Boundary Hardening 

Voce hardening is not considered in the following, but hardening related to the GBs is 
(cf. NLC35G/NLC5G in Table 1). It can be seen in Fig. 2 and Fig. 3 that, depending on 
the chosen rotation angle ip of the central grain, the strain profiles of the GP simulations 
are matching the DDD profiles better, either in the central grain, or in the bounding 
grains. Nevertheless, in the vicinity of the GBs, the plastic strain evolution is captured 
significantly better compared to the Voce hardening approach. For both simulations with 
GB hardening, identical parameters are used (see Table 1). 

5.3. Gradient Model with Gradient Hardening and without Grain Boundary / Voce Hard¬ 
ening 

It is recalled that in the GP model, gradients enter the theoretical framework via the 
defect energy, see Eq. (10). In order to isolate the influence of the defect energy on the 
evolution of the strain gradients from the misorientation and the interaction behavior of 
dislocations across GBs, an additional case is investigated in the following. Therefore, 
the two bounding grains are set to be purely elastic. Only the central grain is elastic- 
plastic and has an ideal (lOO)-orientation with respect to the a;-axis (cf. Table I). All GB 
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Figure 3: NLC5G: Stress-strain curves of DDD and GP simulations (left). Distribution of 
cross-section averaged plastic strain along loading direction of DDD and GP simulations 
(right), obtained at fixed overall plastic strain values (left). Boundary conditions with 
restricted lateral contraction, central grain rotated 5° around loading axis. GB hardening 
used in GP model, no Voce hardening used. 


contributions and the Voce hardening contribution are neglected. In this case, lateral 
contraction of the bounding planes is not restricted (see [53] for the respective BGs). 
As a consequence of the isolation of strain gradients, information about the quality of 
the used defect energy approach in the GP model is obtained. Gompared to the other 
investigated settings, the elastic overall response of the DDD simulations is harder. This 
is due the confinement of plastic activity to the central grain. Gonsequently, for this 
special case, the Young’s modulus used in the GP model has to be slightly adjusted to 
a value of 69.4 GPa. In all other cases, however, the microplasticity-effect (a seemingly 
smaller Young’s modulus) is more pronounced because a small bow out of a favorably 
oriented FR source at the beginning of the simulations is more likely. When GB effects 
and bulk hardening are neglected in the GP model at hand, the defect energy exclusively 
controls the overall rate of work hardening. The plastic defect parameter is fitted to 
a value of Kq = 18 x 10“® N. Fig. 4 shows that the used quadratic form of the defect 
energy leads to an overestimation of plastic strain in the center of the grain. Glose to the 
GBs, the plastic strain is underestimated as are the strain gradients. The overall shapes 
of the plastic strain profiles from both simulation approaches differ more pronounced 
from each other, compared to the cases considered above. 

6. Discussion of Results 

6.1. Grain Boundary Yielding in the Gradient Model 

From the above results, it can be seen that GB yielding is mandatory in the GP model, if 
the plastic strain profiles obtained from DDD simulations with impenetrable GBs should 
be captured. The use of common microhard conditions for the GBs, i.e., an infinitely 
high GB yield strength Sq in the model at hand, is not appropriate. This is due to the 
fact that the DDD profiles show that plastic strain, and thus, plastic slip, is present in 
the immediate vicinity of the GBs. Microhard conditions, however, have been used in 
the GP model of Aifantis et al. [4], where the DDD simulation setup and framework is 
similar to the one in the present work. The evaluation of the plastic strain prohles differs 
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Figure 4: LCOE: Stress-strain curves of DDD and GP simulations (left). Distribution of 
cross-section averaged plastic strain along the loading direction of DDD and GP simula¬ 
tions (right, values are obtained at fixed overall plastic strain values, cf. left). Boundary 
conditions with free lateral contraction, all three grains in (lOO)-orientation. Bounding 
grains behave elastic (i.e., —>■ oo in the GP model of these). GB hardening is neglected 
in the GP model, as are the GB yield strength and the Voce hardening. 


from the present work: there, the DDD strain profiles are obtained by integrating the 
dislocation density tensor with an integration constant, yielding vanishing plastic strain 
at the GB. In the present work, however, the strain profiles are obtained by slicing the 
sample and using the averaging described in Section 3.4, leading to non-zero values of the 
plastic strain in the immediate vicinity of the GBs. Due to the volumetric averaging, 
is evaluated in the slices adjacent to the GBs but not in the GB planes themselves. 

6.2. Grain Boundary Hardening in the Gradient Model 

In order to account for the pile-ups of the DDD results, an additional hardening relation 
is introduced in the GP model. It is associated to the micromorphic field variable ( and, 
consequently, to the equivalent plastic strain 7 eq, at the GB itself. Thus, the GB yield 
strength increases with rising 7 eq. The DDD plastic strain profiles and their evolution 
can not be captured, if, instead of GB hardening, solely a bulk hardening relation of Voce 
type is used in the GP model. This is especially apparent close to the GBs, cf. Fig. 2. 
The DDD results can, furthermore, be captured notably well close to the boundaries 
with this procedure, for the investigated cases at hand. Consequently, this indicates that 
hardening associated to the boundaries as well is necessary in the GP model to account 
for the BCs of the DDD simulations. 

In addition, the use of a finite GB strength Sq , in combination with a sufficiently large GB 
hardening parameter IFh oo, reproduces GB behavior similar to microhard conditions 
(cf. Fig. A.l). However, in this combination, the model behavior allows for some plasticity 
to occur at the GB, at first. Subsequently, by the high GB hardening contribution to 
the GB yield strength, plasticity is prevented from proceeding further. The deviations of 
this approach from ideal microhard conditions (i.e., Sq —>■ oo) are small for the overall 
mechanical responses and for the plastic strain profiles . 
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6.3. Dislocation Interaction across Grain Boundaries 

It is challenging to incorporate dislocation interaction across GBs in GP simulations 
with the same physical details as in DDD simulations. With the used GP model in this 
work, identical plastic strain profiles are obtained for two different crystal orientations 
of the central grain (cf. Fig. 2 and Fig. 3). Contrarily, the DDD simulations show 
significantly different plastic strain profiles due to the orientation dependent dislocation 
interaction across GBs. While no dislocations are allowed to pass the GBs in the DDD 
setting, the stress fields of dislocation pile-ups on adjacent sides of the GBs interact 
nevertheless. This influences the plastic strain gradients close to the GBs, depending 
on the misorientation between the grains. A high misorientation, resulting in a large 
mismatch of the respective stress fields, induces a parabolic distribution of plastic strain 
(cf. Fig. 2). In this case, the interaction between dislocations across the GBs is weaker 
than it is for lower misorientations. Thus, the central grain becomes the main carrier of 
plasticity in the DDD simulations. Regardless of the phenomenological approach to GB 
strengthening, the plastic strain profile in the central grain is reasonably captured by the 
GP model. Gontrarily, a low misorientation leads to a more homogeneous distribution of 
the plastic strain in the grains with sharper “cusps” towards the GBs (cf. Fig. 3). In this 
case, the DDD and GP profiles are in rather good agreement in the two bounding grains. 
Thus, further research regarding the appropriate modeling of the dislocation interaction 
across the GBs in the continuum model is necessary. 

6.4. Evolution of Plastic Strain Gradients 

The gradients of the plastic strain, close to the GBs in the central grain of the DDD 
results, are reasonably captured by the first GP strain profile at Sp = 0.001 (cf. Fig. 2 
and Fig. 3). Furthermore, it can be observed that the GP model with GB hardening 
does capture the evolution of plastic strain values close to the GB, but not the evolution 
of the gradients, there. One could argue that this might stem from the employed GP 
model with 7 eq, (, and its gradient, in general. The investigated special case LGOE, 
however, gives indications that the differences in the gradients between GP and DDD 
results could mainly be caused by the convenient choice of a quadratic defect energy. A 
more general form of the defect energy that could also be formulated with respect to VC 
(e.g., using a power law approach [10]) could give better agreement in the gradients. This 
could also be obtained by applying more physically enriched theories, e.g., [11, 12, 32], 
allowing to consider dislocation densities and a physically more sound defect energy 
definition, respectively. One might further argue that the evaluation of strain gradients 
is dependent on the discretization, especially in the DDD setting. It can, however, be 
shown that the plastic strain values close to the GB encounter only minor differences upon 
further discretization refinement in the DDD simulations. This is due to the fact that 
the discretization resolution is chosen higher than the dislocation spacing in the pile-ups 
(cf. Section 3.4). The GP strain profiles show more pronounced hyperbolic distributions 
than the DDD strain profiles do, independent of the chosen misorientations. The DDD 
strain profiles (and gradients), thus, are obviously sensitive to the misorientation and 
show the role of dislocation interaction across the GBs, a feature, the GP model at 
hand cannot account for. This shortcoming could be removed by using, e.g., one of the 
above-mentioned approaches of [11, 12, 32]. 
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6.5. Assumptions in the Gradient Plasticity Model 

In the GP model, an equivalent plastic strain is taken into account. Its micromorphic 
counterpart ( is assumed to be continuous across the GBs while the gradient is not. 
It is remarked that employing continuity of, e.g., plastic strain, but allowing for jumps in 
its gradient, has shown good agreement with the size effects observed in experiments, in a 
one-dimensional case (cf. [5] , and the discussion in the introduction of the work at hand). 
Furthermore, the equivalent plastic strain of the present work is the (directionally inde¬ 
pendent) sum of the slip contributions of all slip systems. Therefore, it is not a quantity 
from which conclusions on the behavior of individual slips should be drawn, but, rather, 
on an overall response of these. Gonsequently, the choice of continuity for ( across GBs 
is not expected to be in general transferable to statements about the continuity of slips 
on individual slip systems across GBs. These are known to be often discontinuous across 
GBs due to different operating mechanisms during slip transfer (e.g., [38]). Within the 
context of an overall description of plastic slip, it might also be noteworthy that the cal¬ 
culations of an effective plastic strain in the experimental work of [I] lead to a continuous 
appearance of distributions of this quantity across many GBs. Furthermore, the DDD 
simulations in the work at hand show similar values of the plastic strain on the adjacent 
sides of the GBs in a;-direction. As there is a formal connection between 7 eq (and thus 
Q and the components of e^, this finding might indicate that the continuity assumption 
on C could be justifiable for certain cases. 

The internal length scale in the GP model can be deduced from the defect energy 
parameter Kq. The used value of Kq (cf. Section 4.2) yields an internal length of 
I — y/Kc/E « 36 nm, considering a Young’s modulus of 65 GPa. Remarkably, this re¬ 
sult is of the same order of magnitude as the mean dislocation spacing in the pile-ups 
at the GBs of the DDD simulations (cf. Section 3.4). Additionally, it is not expected 
that Kq (and thus, the internal length scale) being a constant should be generally valid. 

7. Conclusions 

For the modeling of pile-ups, as they are observed in discrete dislocation dynamics sim¬ 
ulations with impenetrable grain boundaries, grain boundary yielding is used in the 
gradient plasticity model of this work. The increasingly entangling dislocations at the 
grain boundaries are modeled by using grain boundary hardening.. Thereby, the over¬ 
all mechanical response of discrete dislocation dynamics simulations as well as the local 
response in form of plastic strain profiles is reasonably captured for certain cases. The 
comparison of this approach to a previous model without grain boundary hardening, but 
with an isotropic Voce hardening relation instead, shows inferior agreement. Without 
the grain boundary hardening, the evolution of plastic strain values, observed in the 
immediate vicinity of grain boundaries in discrete dislocation dynamics simulations, can 
not be captured by the continuum model. Microhard grain boundary conditions, com¬ 
monly used in continuum models in the literature, are not sufficient to account for the 
observed behavior at the impenetrable grain boundaries of discrete dislocation dynamics 
simulations. Gontrary to these common conditions, rather a finite grain boundary yield 
strength needs to be used. The plastic strain profiles obtained from the continuum model 
are in agreement with the discrete dislocation dynamics results in the immediate vicinity 
of the grain boundary. It could also be argued that physically more advanced theories (by 
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incorporating higher-order terms) might account for the sharp gradients close to grain 
boundaries - in spite of using microhard conditions. The present GP model, however, al¬ 
lows for extended fully threedimensional simulations in manageable computational times, 
a benefit often lost when using physically enriched theories with substantially increased 
degrees of freedom. Differences between the results of the continuum and the discrete 
modeling approach indicate that additional research is necessary related to accounting 
for crystal orientation dependent dislocation interaction mechanisms across grain bound¬ 
aries in the continuum model. This also includes the possible recalibration of the model 
parameters for different orientations and loading conditions as well as further case stud¬ 
ies. Additionally, the question of an enhanced defect energy approach is raised, as well 
as a suitable combination of hardening approaches considering both, geometrically nec¬ 
essary dislocations and statistically stored dislocations. Further investigations could also 
consider elastic anisotropy and different boundary conditions. 


Appendix A. Parameter Study of Grain Boundary Hardening 
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Figure A.l: Parameter study of grain boundary hardening parameter Ah. Case NLC35G, 
only Ah is varied as shown above. All other parameters can be found in Table 1 and 
Section 4.2, respectively. 
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